{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "8.625716949943694\n",
      "9.270944409698814\n",
      "9.96443666613547\n",
      "10.709804059395822\n",
      "11.510926993039504\n",
      "12.371976135533515\n",
      "13.297434132869398\n",
      "14.292118945343816\n",
      "15.361208929995087\n",
      "16.51026979927539\n",
      "17.745283596305597\n",
      "真实值: [ 7.19  8.26  9.5   9.98 10.48 11.   18.13 16.46 14.79 15.66 16.3  16.66]\n",
      "预测值: [19.1 20.5 22.  23.7 25.5 27.4 29.4 31.6 34.  36.5 39.2 42.2]\n",
      "{'a': {'value': -0.07213716836353937, 'desc': '发展系数'}, 'b': {'value': 8.421908304326227, 'desc': '灰色作用量'}, 'predict': {'value': array([19.07267984, 20.49936899, 22.03277842, 23.68089112, 25.45228721,\n",
      "       27.35618861, 29.4025071 , 31.60189587, 33.96580499, 36.50654105,\n",
      "       39.23733118, 42.17239196]), 'desc': '往后预测12个的序列'}, 'C': {'value': 0.25974677344500624, 'desc': '后验差比<=0.35，模型精度等级为好'}}\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAydUlEQVR4nO3dd3hUZdrH8e9NUaTqIiAgGNZFEelEOiogRURcRJYi1hWw42ID9V2aq+jaV0VQWFARXFdQUUQQRJoCQUFROkYXQao0kZY87x/PJAxhAoFkciaT3+e65srMmXPO3AmH3DlPuR9zziEiIpJRgaADEBGR2KQEISIiESlBiIhIREoQIiISkRKEiIhEVCjoAHLSmWee6RISEoIOQ0Qkz1i8ePFW51yZSO/FVYJISEggKSkp6DBERPIMM/sxs/fUxCQiIhEpQYiISERKECIiElFc9UFEcvDgQdavX8++ffuCDkVOQpEiRTj77LMpXLhw0KGI5DtxnyDWr19PiRIlSEhIwMyCDkdOgHOObdu2sX79eqpUqRJ0OCL5Ttw3Me3bt4/SpUsrOeRBZkbp0qV19ydyHIMGRee8cZ8gACWHPEz/diLHN3hwdM6bLxKEiEi8+vLL6J1bCSLKduzYwcsvvxx0GOnGjBnDnXfeGXQYIpJNgwaBGTRu7F+b+UdONjcpQUSycSNccgn88ku2T3WsBJGSkpLt84tI/nTFFVCyJJx7rn/tnH8oQUTb0KEwdy4MGZLtU/Xv35+1a9dSp04d7r//fmbNmkWLFi3o0aMHNWvWJDk5mRo1aqTv/9RTTzEo9C+8du1a2rVrR/369WnevDkrVqw44typqakkJCSwY8eO9G1/+tOf2LRpE5MnT6Zhw4bUrVuXyy67jE2bNh0V24033sh///vf9NfFixdPf/7Pf/6Tiy66iFq1ajFw4MBs/xxEJOcsXgxt2kDp0vDZZ9H7nLgf5nqEe+6BJUsyf3/OHEhNPfx6+HD/KFAAmjePfEydOvDcc5mectiwYSxbtowloc+dNWsWCxcuZNmyZVSpUoXk5ORMj+3duzevvPIKVatWZcGCBdx+++3MnDkz/f0CBQpw1VVXMWnSJG666SYWLFhAQkIC5cqVo1mzZnz55ZeYGa+99hpPPvkkTz/9dObfe5hp06axevVqFi5ciHOOjh07Mnv2bC6++OIsHS8i0fPVV9C6NZx+uk8OlSpBtP6Gy18J4ngaNIB162DrVp8oChSAM888fA+XYx/T4Ljj+vfs2cP8+fPp0qVL+rb9+/cftV/Xrl0ZMmQIN910ExMmTKBr166An//RtWtXNm7cyIEDB05oHsG0adOYNm0adevWTY9l9erVShAiAVuyBC67DEqU8MnhnHP89mgNc41agjCzSsDrwFlAKjDSOfe8mQ0CegFbQrs+5JybEuH4dsDzQEHgNefcsGwHdYy/9NPddhuMHAlFisCBA9C5M+RwJ3OxYsXSnxcqVIjUsLuWtDH/qampnH766el3Hplp3Lgxa9asYcuWLbz33ns88sgjANx1113069ePjh07MmvWrPRmq3Dhn+2c48CBA+nPBwwYQJ8+fbLzbYpIDlq6FFq1guLFYdYsyI2VDaLZB3EIuNc5dwHQCLjDzKqH3nvWOVcn9IiUHAoCLwGXA9WB7mHHRtemTXDrrX7s2K23ZrujukSJEuzevTvT98uVK8fmzZvZtm0b+/fv58MPPwSgZMmSVKlShXfeeQfwv7SXLl161PFmRqdOnejXrx8XXHABpUuXBmDnzp1UrFgRgLFjx0b87ISEBBYvXgzA+++/z8GDBwFo27Yto0ePZs+ePQD8/PPPbN68+WS+fRHJAd9+65ND0aL+ziG3CgtE7Q7CObcR2Bh6vtvMlgMVs3h4A2CNc24dgJlNAK4Cvo9GrEeYOPHw85deyvbpSpcuTdOmTalRowaXX345V1xxxRHvFy5cmL///e80bNiQKlWqUK1atfT3xo0bx2233cajjz7KwYMH6datG7Vr1z7qM7p27cpFF13EmDFj0rcNGjSILl26ULFiRRo1asQPP/xw1HG9evXiqquuokGDBrRq1Sr9zqZNmzYsX76cxqHxc8WLF+fNN9+kbNmy2f55iMiJ+e47nxyKFPHJIYdbvI/JnHPR/xCzBGA2UAPoB9wI7AKS8HcZv2bY/xqgnXPultDr64CGzrmjBvCbWW+gN0DlypXr//jjkWtfLF++nAsuuCCHvyPJTfo3lPzq+++hRQsoWBA+/xyqVs35zzCzxc65xEjvRX2Yq5kVB94F7nHO7QKGA+cCdfB3GJGG1kSqrxAxkznnRjrnEp1ziWXKRFw1T0Qkz1mxAlq29GNlPvssOsnheKKaIMysMD45jHPOTQRwzm1yzqU451KBV/HNSRmtByqFvT4b2BDNWEVEYsXKlf7OAXxyOP/8YOKIWoIwX2VtFLDcOfdM2PbyYbt1ApZFOHwRUNXMqpjZKUA34INoxSoiEitWrfLJITUVZs6EsG7JXBfNeRBNgeuAb81sSWjbQ/gRSXXwTUbJQB8AM6uAH87a3jl3yMzuBD7BD3Md7Zz7LoqxiogEbs0anxwOHfJ3DtVzZ+xmpqI5imkukfsSjhrWGtp/A9A+7PWUzPYVEYk3a9f65HDggL9zuPDCoCPSTGoRkcCtW+eTw++/w4wZULNm0BF5KtaXx8yaNYsOHToA8MEHHzBsWOYTzE+21PigQYN46qmnTjrGnD6PSDxLTvbJ4bff4NNPIcJUp8AoQWQiWrVNMnMypb87duxI//79M30/1taiEJHDBg2CH3/0yWHXLpg+3df+jCVKEJnIqSX8kpOTqVatGjfccAO1atXimmuuYe/evYAvdTFkyBCaNWvGO++8w7Rp02jcuDH16tWjS5cu6aUupk6dSrVq1WjWrBkTw2Z6hy/+s2nTJjp16kTt2rWpXbs28+fPP6rUOGRexvsf//gH559/PpdddhkrV6486vvYuXMnCQkJ6bWb9u7dS6VKlTh48CCvvvoqF110EbVr16Zz587p31+4Sy+9lKSkJAC2bt1KQqiQTEpKCvfff396TCNGjMjuj1wkTxg82CeHHTv8nUO9ekFHdLR81QdxvGrfGV166fH3OU61bwBWrlzJqFGjaNq0KTfffDMvv/wy9913HwBFihRh7ty5bN26lauvvppPP/2UYsWK8cQTT/DMM8/wwAMP0KtXL2bOnMmf/vSn9GqtGd19991ccsklTJo0iZSUFPbs2XNUqfHMyngXK1aMCRMm8PXXX3Po0CHq1atH/fr1jzh/qVKlqF27Np9//jktWrRg8uTJtG3blsKFC3P11VfTq1cvAB555BFGjRrFXXfddfwfHjBq1ChKlSrFokWL2L9/P02bNqVNmzYnVH1WJK9Zv95/3b7d3zlk+O8WM3QHESY52U9n//xz/zrt+TGWbMiSSpUq0bRpUwB69uzJ3Llz099L+4X/5Zdf8v3339O0aVPq1KnD2LFj+fHHH1mxYgVVqlShatWqmBk9e/aM+BkzZ87ktttuA6BgwYKUKlXqqH3Cy3jXq1ePFStWsHr1aubMmUOnTp0oWrQoJUuWpGPHjhE/o2vXrrz99tsAR5QWX7ZsGc2bN6dmzZqMGzeO777L+ojkadOm8frrr1OnTh0aNmzItm3bWL16dZaPF8lL0pYJrRSaBrxzp19lILebtLMqX91BZKXadxozv3xfTvBzBiO/TiuQ55yjdevWjB8//oh9lyxZctTxJyuzMt7PPfdclj6jY8eODBgwgO3bt7N48WJatmwJ+JXp3nvvPWrXrs2YMWOYNWvWUceGlxZPK2meFtO//vUv2rZtm43vTCRv6N4dXnsN9u6FX3/Nud8x0aI7iFzw008/8cUXXwAwfvx4mjVrdtQ+jRo1Yt68eaxZswbwbfyrVq2iWrVq/PDDD6xduzb9+EhatWrF8OHDAd+uv2vXrqNKjWdWxvviiy9m0qRJ/P777+zevZvJkydH/IzixYvToEED+vbtS4cOHShYsCAAu3fvpnz58hw8eJBx48ZFPDa8tHj4Mqdt27Zl+PDh6aXGV61axW+//RbxHCJ52TffwMUXw8GDfj2HvEAJIhM5uYTfBRdcwNixY6lVqxbbt29PbwoKV6ZMGcaMGUP37t2pVasWjRo1YsWKFRQpUoSRI0dyxRVX0KxZM85JW0Iqg+eff57PPvuMmjVrUr9+fb777rsjSo3ff//9tGnThh49etC4cWNq1qzJNddcw+7du6lXrx5du3alTp06dO7cmeaZLa+Kb2Z68803j+gLGTp0KA0bNqR169ZHlCsPd9999zF8+HCaNGnC1q1b07ffcsstVK9enXr16lGjRg369OnDoUOHsvqjFckTFi3yfZqFC8Ps2VCrVvSWCc1JuVLuO7ckJia6tJEyaYIuFZ2cnEyHDh1YtixSySnJiqD/DUWyY/Zs6NDBr148Y0buLfaTVYGW+xYRya+mTYN27aBiRZgzJ/aSw/EoQURZQkKC7h5E8qH334crr4TzzvOjIStmdT3NGJIvEkQ8NaPlN/q3k7xo/Hjo3Bnq1vVVWfPqar1xnyCKFCnCtm3b9IsmD3LOsW3bNooUKRJ0KCJZNno0XHstNG3qJ8GdcUbQEZ28uJ8HcfbZZ7N+/Xq2bNkSdChyEooUKcLZZ58ddBgiWfLCC9C3L7RtCxMnQtGiQUeUPXGfIAoXLqyyDSISdY8/Dg89BJ06+SamU08NOqLsi/smJhGRaHIOHn7YJ4drr4X//Cc+kgNEd03qSmb2mZktN7PvzKxvaPs/zWyFmX1jZpPM7PRMjk82s2/NbImZJUXaR0QkSM7B3/4Gjz0GvXrB2LFQKI7aZaJ5B3EIuNc5dwHQCLjDzKoD04EazrlawCpgwDHO0cI5VyezSRwiIkFJSYHeveH5532/w4gREKo+EzeiliCccxudc1+Fnu8GlgMVnXPTnHNptRS+BNQDKSJ5ysGDcP31vvDeI4/As8/6Ap/xJlf6IMwsAagLLMjw1s3Ax5kc5oBpZrbYzHof49y9zSzJzJI0UklEom3/fvjLX+Ctt3zH9NCh8ZkcIBdGMZlZceBd4B7n3K6w7Q/jm6Eil/+Eps65DWZWFphuZiucc7Mz7uScGwmMBF+LKce/ARGRkIcfhsWL4ZNP/JDWLK6LlWdFNUGYWWF8chjnnJsYtv0GoAPQymUyg805tyH0dbOZTQIaAEclCBGR3LBrl++MNoNRo+Dmm4OOKPqiOYrJgFHAcufcM2Hb2wEPAh2dc0cvXuz3KWZmJdKeA20AFTQSkUBs3gyh9bF46638kRwgun0QTYHrgJahoapLzKw98CJQAt9stMTMXgEwswpmNiV0bDlgrpktBRYCHznnpkYxVhGRiO6+G8qV801L4FeFM4vdZUJzUtSamJxzc4FIXTdTImxLa1JqH3q+DqgdrdhERLJi8WJ4+20oXRo+/BAaN479ZUJzkmZSi4hEMH26XwXutNNg3jxo1CjoiHKfEoSISAbjxkH79nDuuTB/Ppx/vt+eF5YJzUlKECIiYZ5+Gnr2hObN/UI/FSocfi8/9DuEU4IQEQFSU6FfP7jvPj8R7uOPoVSpoKMKlhKEiOR7+/f7u4Znn/WjluKlXHd2xVHdQRGRE7drF1x9NcyYAU88AfffH7+lM06UEoSI5Fu//AKXXw7LlvlS3ddfH3REsUUJQkTypVWroF07P0t68mT/XI6kBCEi+c7ChXDFFb4p6bPP4KKLgo4oNqmTWkTylY8/hhYtoEQJPwFOySFzShAikm+MHQtXXuknvs2fD1WrBh1RbFOCEJG45xwMGwY33ujvHj7/HM46K+ioYp/6IEQkrv3977BjB/zrX9CjB/z733DKKUFHlTcoQYhI3Nq3zy8JCnDvvfDkk1BA7SZZph+ViMSlrVuhTRv//Kmn/EPJ4cToxyUiceeOO6BMGZgzx7++7778s8hPTlITk4jElWnTfLnusmXh/ffz3yI/OUl3ECISN15+2a/jULmynwyXHxf5yUlRSxBmVsnMPjOz5Wb2nZn1DW3/g5lNN7PVoa9nZHJ8OzNbaWZrzKx/tOIUkbzv0CFfhfWOO3xtpXnz4Jxz/Hv5bZGfnBTNO4hDwL3OuQuARsAdZlYd6A/McM5VBWaEXh/BzAoCLwGXA9WB7qFjRUSOsHMndOjgh7Heey+8956fJZ1G/Q4nL2oJwjm30Tn3Vej5bmA5UBG4Chgb2m0s8OcIhzcA1jjn1jnnDgATQseJiKT74Qdo0sSX6h450o9UKlgw6KjiR650UptZAlAXWACUc85tBJ9EzKxshEMqAv8Le70eaJjJuXsDvQEqV66cg1GLSCybOxc6dYKUFN8x3aJF0BHFn6h3UptZceBd4B7n3K6sHhZhW8RxCM65kc65ROdcYpkyZU42TBHJQ954A1q1gjPOgC+/VHKIlqgmCDMrjE8O45xzE0ObN5lZ+dD75YHNEQ5dD1QKe302sCGasYpI7EtNhYcf9gv7NG3qk8N55wUdVfyK5igmA0YBy51zz4S99QFwQ+j5DcD7EQ5fBFQ1sypmdgrQLXSciORTe/fCX/4Cjz0GvXrBJ5/AH/4QdFTxLZp3EE2B64CWZrYk9GgPDANam9lqoHXoNWZWwcymADjnDgF3Ap/gO7f/45z7LoqxikgM27ABLr4YJk6EZ56BESOgcOGgo4p/Ueukds7NJXJfAkCrCPtvANqHvZ4CTIlOdCKSV3z1FXTs6Iezvv++X89BcodmUotIzJo0CZo390NX581TcshtqsUkIjFn4EAoWhT694eGDf3kNy3wk/uUIEQkpuzfD0OG+OfdusHo0XDaacHGlF+piUlEYkZyMlx6qX8+aBC89ZaSQ5CUIEQkJnTpAlWq+LkN4BNEgQKqpRQkNTGJSKD27oW+feG///Xlud96C/74R63hEAt0ByEigfnmG0hMhFGjYMAAmD3b30VIbNAdhIjkOudg+HDo18/XU5o+3ddWSqM1HGKD7iBEJFdt3w6dO/vFfVq2hKVLj0wOoH6HWKEEISK5Zs4cqF0bPvwQnn7afy0bqeC/xAQlCBGJupQUGDzYD2EtUgS++MI3LxXQb6CYpj4IEYmq9evh2mt9B/R118FLLx25JKjELiUIEYma99+Hm2/2s6Nff90nCMk7dIMnIjlu3z646y74858hIQG+/lrJIS9SghCRHLV8uS+w9+KLvp9h/nyoWjXoqORkqIlJRHLEwIFQuTLcfTcUKwYffQTt2x//OIldShAikm07dx6uwNqqFbzxBpQvH2xMkn1qYhKRbFmwAOrW9c8fe8yvFa3kEB+iliDMbLSZbTazZWHb3g5bnzrZzJZkcmyymX0b2i8pWjGKyMlLTYXLLvMF9n74wW976CEoVEgzoeNFNJuYxgAvAq+nbXDOdU17bmZPAzuPcXwL59zWqEUnIiftl1/8qKQZM+Avf4ERI3xNJVVgjS9RSxDOudlmlhDpPTMz4C9Ay2h9vohEx9SpcP31sGcPvPoq/PWvYBZ0VBINQfVBNAc2OedWZ/K+A6aZ2WIz632sE5lZbzNLMrOkLVu25HigIuIdOAD33QeXX+7Xh05KgltuOZwcVIE1/hw3QZhZOTMbZWYfh15XN7O/ZvNzuwPjj/F+U+dcPeBy4A4zuzizHZ1zI51zic65xDJlymQzLBGJZM0aaNLEF9i7/XbfMV29+pH7qN8h/mTlDmIM8AlQIfR6FXDPyX6gmRUCrgbezmwf59yG0NfNwCSgwcl+nohkz5tv+lFK69bBpEm+lpLWic4fspIgznTO/QdIBXDOHQJSsvGZlwErnHPrI71pZsXMrETac6ANsCzSviISPXv2wA03+M7ounX9ug1//nPQUUluykqC+M3MSuP7BTCzRhx79BGh/cYDXwDnm9n6sGapbmRoXjKzCmY2JfSyHDDXzJYCC4GPnHNTs/TdiEiO+OorqFfP3z0MHAgzZ0KlSkFHJbktK6OY+gEfAOea2TygDHDN8Q5yznXPZPuNEbZtANqHnq8DamchLhHJYc7B88/DAw/4hXxmzoRLLgk6KgnKce8gnHNfAZcATYA+wIXOuW+iHZiI5K4tW6BDB/jb33wNpaVLTyA5bNzod/7ll6jGKLkrK6OYrgd6APWBekD30DYRyePSRh7NnOmXAp0xw1dhnTQJSpc+gRMNHQpz5x4uyCRxIStNTBeFPS8CtAK+ImyGtIjkTYMHw8GD8PjjcP758PHHPlEc1/btsHixv9U4dOjw9uHD/aNIEfj996jFLbnjuAnCOXdX+GszKwW8EbWIRCRXJCf7r4895mdDP/+8L9N9lB07fK91UpJ/LF7sx7ymKV7cJ4OUFChaFDp1gqeeyoXvQKLtZEpt7AW0/IdIHjVw4NEtQaNGwdlnw6B+u/zyb2nJICnJz5JLk5AAiYnQu7f/Wq+er9A3cqS/a9i3D0qW9FOtJc87boIws8mEhrji+yyqA/+JZlAiEh3r1/tZ0AAtyq/gs43VcE2b+V/8E5JgyKrDFfcqV/ZJ4Kab/Nf69SN3TGzaBLfe6pPGyJG+w1rigrnjlF80s/BxDIeAHzOb5Ba0xMREl5Sk6uAiGTkHr78Offs6Du7cy5M8wG0MpyCpOELFlAoU8J0S9ev7R9mywQYtucLMFjvnEiO9l5U+iM9zPiQRyS0bNzj6dN/F5NmlaH7KAv7NtZxb6CdwBRiYMsjXzbj6at9voKYhCZPpMFcz221muyI8dpvZrtwMUkROnPvxJ97q9gEXVt7F9Nmn8GyBe5nV7gnOfecJuPFGcI5BRZ6A/fvVbyARZXoH4ZwrkZuBiEgO2LED/vtfNo/+kNu+uI6JdKZRiWWM6fct59/10OE+hLfeUr+BHNdx+yDSdzQri58HAYBz7qdoBXWy1AchcW/jRujWDd5++/Bf/Pv3+wkMb74Jkyfz3wNXclvBkeyiJEPv3cm9j5WmYMFgw5bYla0+CDPrCDyNL/e9GTgHWA5cmJNBikgWpM1YHjwYevSAcePgP/+BX39l25nnc0fCF7y9qh6JdR1jxxrVq5/IdGiRI2VlHsRQoBHwqXOurpm1wC/4IyK55bTT/ByDNK+84h8A117L++f2o8+Iumz/wXj0UXjwQaNQNFecl3whK+W+DzrntgEFzKyAc+4zoE50wxKRI0ydCuecc/h1gQLQpAm/Ll7H9QXe5M9D6lG+vJGUBA8/jJKD5IisJIgdZlYcmAOMM7Pn8fMhRCTali6FLl2gRQv4+WcABhUcCsCU03tQ48oqjB/vZ0cvWAC1agUZrMSbrPydMRs4HegL9ARKASrZKBJNSUm+v+GDD/wQ1IcegiVL4JxzGPzyI/yvWktGT2lCjRowebKveCGS07KSIAy/JvV2YALwdqjJSURy2vz5PjFMnQqnn+7rcd99N5xxBuAHK/EyjFnVhAED/J3DqacGGbDEs6wsGDTYOXchcAd+JNPnZvbp8Y4zs9FmttnMloVtG2RmP5vZktCjfSbHtjOzlWa2xsz6n8D3I5I3ff45tGoFTZv6u4fHHoMff/QZ4IwzuP12MPPVtQFSU32J7scfDzZsiW9Z6YNIsxn4BdgGZKVIyxigXYTtzzrn6oQeUzK+aWYFgZeAy/GFAbubWfUTiFMkb3AOpk+Hiy+GSy+F777z5S6Sk2HAAChZMr0O3ogR/obi6acPH+rc4QV/RKIhKyvK3WZms4AZwJlAL+fccbvCnHOz8c1SJ6oBsMY5t845dwDfrHXVSZxHJDY5Bx99BI0bQ5s2fm2FF16AH36Ae++FYsXYtw+GDYOqVX0p7jvv9FW3+/ULOnjJT7JyB3EOcI9z7kLn3EDn3PfZ/Mw7zeybUBPUGRHerwj8L+z1+tA2kbwpbb3mDRv8Wp6JiX7x519+8XMZ1q6Fu+6C007DORg/HqpV8zcRLVv6G4vnnz9cJWPgwGC/Hck/stIH0d85tySHPm84cC5+HsVG/AztjCxSGJmd0Mx6m1mSmSVt2bIlR4IUyVGDB8OcOVCzpq+aumsXjB4Nq1dDnz7pvczz5/ubih49fJ/0zJnw3ntw3nlHnk7NSpJbcnU6jXNuU9pzM3sV+DDCbuuBSmGvzwY2HOOcI4GR4Gsx5UykIjkg4+zn7aEW1/Xr/SI8IevWQf/+8M47UKEC/PvfcN11qH6SBO5EOqmzzczKh73sBCyLsNsioKqZVTGzU4BuwAe5EZ9Ijtm6Fa4KdZ1Z6Kb4tNPg2mt9XwO+8OoDD8AFF/guiUGDYNUqX4lbyUFiQdTuIMxsPHApcKaZrQcGApeaWR18k1Ey0Ce0bwXgNedce+fcITO7Ez/3oiAw2jn3XbTiFMlRqam+V7l/f9+UVKcOfPONb0YKrbtwsPRZjHzJJ4Rt2+CGG+DRR6GietokxkQtQTjnIhX0G5XJvhuA9mGvpwBHDYEViWlLl8Jtt8EXX0Dz5vDyy/D3v0OTJgxK+T8GFhjKlCUVuK8WrFjhq2c8/TTUrRt04CKRqaSXSHbt3u2HFr3wgu9dHjMGrr/eNy1NnAjAYIO5rV5ixhe+0/n99+HKKw+3PonEolztgxCJK875nuVq1eC55+CWW2DlSt9mFPrNv3Gj3wzw9dc+hyxbBh07KjlI7NMdhMjJWLPGz1775BPfz/Duu9CoUfrbe/fCFVfArFmHD9m+3ZdV2rZNQ1Ulb9AdhMiJ2LfPz2uoUcNPXHjuOVi0KD05pKbCG2/4ZqRZs6BzZ59LQOUxJO/RHYRIVk2bBnfc4X/jd+0KzzzjJy6EfP65r5SxeLGfLD1+vO+rFsmrdAchEklaeYxffvEL9XTtCm3b+o6DadNgwoT05LB6tZ8gfemlsGkTvPmmX7wnPDmoPIbkRbqDEIlk6FCYOxeuucYPXz140DctPfAAFCkC+D6FoUPhpZf8NIdHH4W//Q2KFj36dGpWkrxICUIkXMbyGPPm+a+nnurnNAAHDsDw4T5f7NwJf/0rDBkCZ50VQLwiUaQmJpFw33/va1+kOeUUXz0vORnnfPG8Cy+Ee+7x/QxLlsDIkUoOEp+UIETSLFrkx6YuX+5fFykChw5BqVJ8teEsWrSATp2gcGGYMsWPcK1ZM9iQRaJJCUJk/3545BFfa3vXLmjSBG6/nUE3/MD6nv254cMuJCb6m4uXX/allS6/XBPdJP6pD0Lyt6+/9jOfv/3Wf33uOTj9dPbsgcEl4MnT/kFKiu+bHjAASpUKOmCR3KM7CMmf0kYlNWgAW7b44khjxpBS4nRGjz68SE/Hjr6w3rBhSg6S/yhBSP7z7bfQsKEfe/qXv6QXR5oxw5fc/utf/TQIgLffhj/+UcNUJX9SgpD849AhePxxqF/fr+r27rswbhwrtpTmyivhsst8v/SECb5kBqg8huRv6oOQ/GH5cr9U28KFfvLbyy+z1cow+C4/p6FoUd+M1Ldv+jw4kXxPdxAS31JS4Kmn/Ko8a9bAhAnsf/Mdnhpbhj/9ySeH3r39Ww8+eGRyUHkMye90ByHxa/Vqf9cwfz5cdRVu+Cu8O+8sHrjALwvdvj38859QvXrkw9WsJPld1O4gzGy0mW02s2Vh2/5pZivM7Bszm2Rmp2dybLKZfWtmS8wsKVoxSpxJK7C3YQM8/zzUru0nL7zxBgsHTKJ5l7Po0gWKF/f19j76KPPkICLRbWIaA7TLsG06UMM5VwtYBQw4xvEtnHN1nHOJUYpP4s3QoTBnjq+Bcc890KIFP01bwbUf96RhI2PNGnj1VT/1oXXroIMViX1Ra2Jyzs02s4QM26aFvfwSuCZany/5SFiBvUEMZNDGweyiBMM+uYRnZ5YD4OGHfR9DiRJBBiqStwTZSX0z8HEm7zlgmpktNrPexzqJmfU2syQzS9qyZUuOByl5wKpV6TPbBjOIkYVup2qR//F4ygNcc41/+9FHlRxETlQgndRm9jBwCBiXyS5NnXMbzKwsMN3MVjjnZkfa0Tk3EhgJkJiY6KISsMSunTuhVy/cqlVMDbVo9jn0Es3LrObD90tx0UUBxyeSh+X6HYSZ3QB0AK51zkX8he6c2xD6uhmYBDTIvQglz/jhB1zjJlz3SU8K4GgfdkM6Z2NVPvoowNhE4kCuJggzawc8CHR0zu3NZJ9iZlYi7TnQBlgWaV/Jv9ycucyocy/NV43iTXpSqRK88kroPc1+FskR0RzmOh74AjjfzNab2V+BF4ES+GajJWb2SmjfCmY2JXRoOWCumS0FFgIfOeemRitOyXtmPTSNSy9J5bJdE0k+sz4vv+ynPPTpE3RkIvElmqOYukfYPCqTfTcA7UPP1wG1oxWX5F2zZ6Uy8IZkZv3UhgqnbOVfQ37jlr7FNPtZJEpUakNi3rx5cFmLFC5pUYAVP53Gc40nsGZLKe58sNhRdZPUrCSSc5QgJGZ9+SW0bQvNmsG3c3fwDP1YO+y/9J3XldNKFg46PJG4pwQhMWfhQr+kZ+PG8PWig/yz5FDWnVqdv33QgqIP3qW1PkVyiRKEBC6tWWjxYujQwa/ls2gRDOu5jHX7KnLf6a9RbP50uPLKQOMUyW9UzVUCN3gwfPUVTJ4MZ5wB/3jUcVfKc5QY2M9ni/feg7POCjpMkXxHdxASmEOHoHNn/3zOHF9rL3nlfh5afZNPDt26wWefKTmIBEQJQgIxaBAULgwTJ/rXO3bA//0fPJP4Fowd63d46y1fiE9EAqEEIYEYNCg043nDRgDc57NxfzyXQZtug/Hj/YQGdUaLBEoJQoI1dKj/2ro1/PYbzJrlm5ZEJHDqpJZghK3hMJCycOAAbNoELVrA778HHJyIgO4gJDft3w/Tp/vV3sqXT988iME+YVx7rV8sWkRigu4gJLp+/hmmTPGP6dN9M9Kpp0LLllChAsyf71/v3w8lS2rEkkgMUYKQnJWSAgsWwEcf+aSwZInfXrkyXH89tG/vk0PRonD11XDbbdC7N4wcCRs3Bhq6iBzJMlmzJ09KTEx0SUlJQYcR/zZu9B3Jb7/t/+Lfvh2mTvUJYepU2LYNChaEpk19QrjiCrjwQo1KEolBZrbYOZcY6T3dQciJGzLEz2z785+hUCH44gtITYUzzzycENq08dOiRSTPUoKQrPn1VyhXDg4ePLxtwQL/tVAh35eQmOjvHEQkLmgUk0SWmuor5g0d6puKzjzTJ4fChQ8ngbSRR//7n6+ZpOQgEleUIOSwTZvgjTf8L/1y5aBBAz+j+eBBePhhv3LPTTf5KdBFimjkkUici+aa1KPNbLOZLQvb9gczm25mq0NfIzZSm1k7M1tpZmvMrH+0YsxXNm6ESy6BX345vO3gQZg9Gx56COrV87/or78ePv3UL8gwbpxPGgsX+n6HJk1gyxa49Va/ms+ttx55PhGJK1EbxWRmFwN7gNedczVC254EtjvnhoV+8Z/hnHsww3EFgVVAa2A9sAjo7pz7/nifqVFMx3D77TBiBPToAc2b+9FGn34Ku3f7PoQmTaBdO/+oXRsK6OZSJD8IZBSTc262mSVk2HwVcGno+VhgFvBghn0aAGucc+sAzGxC6LjjJgiJIK0pKM2bb/oHQJ8+fk3Pli2hVKlg4hORmJXbo5jKOec2AjjnNppZ2Qj7VAT+F/Z6PdAwsxOaWW+gN0DlypVzMNQ8KiUFkpL8rOXp0/3rcIUL+6QwcuQR5S5ERDKKxWGukWZTZdoO5pwbCYwE38QUraBi2tq1hxPCzJl+cQXw/Qr33gvffecnsZ1yii+KV6mSkoOIHFduJ4hNZlY+dPdQHtgcYZ/1QKWw12cDG3IluliTccZymu3bfSJISwppBe4qV/ZLtLVuDa1a+aGp4Eta3HqrSlqIyAnJ7QTxAXADMCz09f0I+ywCqppZFeBnoBvQI9cijCVDh8LcuX6oaffuhxNCUpIfalqypC+P3a+fTwrnnRe5nEXasm0AL72Ue/GLSJ4WzVFM4/Ed0mcCm4CBwHvAf4DKwE9AF+fcdjOrALzmnGsfOrY98BxQEBjtnPtHVj4zbkYxha2VcJSmTX0yaN3az1MoFIuthCKSVxxrFJOK9cWixYt9PaNNm/zrQoV8Ynj1VahaNdjYRCSuHCtBaLB7rJkwAS67zFdENfPDVFNToXp1JQcRyVVqn4gVv/4Kd9wB48dDo0ZQooRPCOpYFpGAKEHEgk8/hRtv9E1Kjz4KDz54ZN+COpZFJABqYgrS77/79Zlbt4bixf26Cg8/rI5nEYkJ+k0UlK+/hp494fvv4a67YNgwvwyniEiM0B1EbktJgccf9+sn7NgBn3wCL7yg5CAiMUd3ELlp3TpfTnvePOjSBV55Bf7wh6CjEhGJSHcQucE5GDXKl9FetsxXU337bSUHEYlpShDRtnkzdOoEt9wCF10E33zjV2yLVBJDRCSGKEHktPCV2yZPhpo14eOP4emn/XBWlSQXkTxCfRA5La3AXsuWsHy5b1aaMQNq1Ag6MhGRE6IEkVMyFthbvtx/XblSyUFE8iQ1MeWUpCQ499zDr0891fc1pK3VICKSxyhBZFdKiq+VdPHFfmU38AX2Dh706zWEL/QjIpKHqIkpOxYt8gX2Fi3yCaJwYTj/fBXYE5G4oARxMrZtg4ce8uszlCvn5zX06HHk0FUV2BORPE5NTCciNdUnhfPO8xPf7rnHd0JrXoOIxKFcTxBmdr6ZLQl77DKzezLsc6mZ7Qzb5++5HedRkpKgcWPffHThhb7Y3jPP+H4GEZE4lOtNTM65lUAdADMrCPwMTIqw6xznXIdcDC2y7dt9Ce4RI6BsWXjjDd0xiEi+EHQTUytgrXPux4DjOFpqKrz2mm9OevVV6NvXNyf17KnkICL5QtAJohswPpP3GpvZUjP72MwuzM2gWLwYmjSBXr38WtBffw3PPgulSuVqGCIiQQosQZjZKUBH4J0Ib38FnOOcqw38C3jvGOfpbWZJZpa0ZcuWkwsmrX7S8uVw++2+qF5ysm9O+vxzX09JRCSfCfIO4nLgK+fcpoxvOOd2Oef2hJ5PAQqb2ZmRTuKcG+mcS3TOJZYpU+bkIhkyBObMgbp1fV/D3XerOUlE8r0gE0R3MmleMrOzzPxvZjNrgI9zW45HcNppPgG88opfs2H/ft/3MGKEmpNEJN8LJEGYWVGgNTAxbNutZnZr6OU1wDIzWwq8AHRzzrkcD2TdOujeHQoW9K+LFlX9JBGRkEBmUjvn9gKlM2x7Jez5i8CLUQ+kfHl/p+Ccr5+0b5/qJ4mIhAQ9iil4mzbBrbfCl1/6r7/8EnREIiIxQbWYJk48/Fz1k0RE0ukOQkREIlKCEBGRiJQgREQkIiUIERGJSAlCREQiUoIQEZGILBoTlINiZluAky0dfiawNQfDySmK68QorhOjuE5MPMZ1jnMuYiG7uEoQ2WFmSc65xKDjyEhxnRjFdWIU14nJb3GpiUlERCJSghARkYiUIA4bGXQAmVBcJ0ZxnRjFdWLyVVzqgxARkYh0ByEiIhEpQYiISERxmyDMrJKZfWZmy83sOzPrG9r+BzObbmarQ1/PCDtmgJmtMbOVZtY2bHt3M/vWzL4xs6mZrY8djbjMrHRo/z1m9mKGc9UPxbXGzF5IW6Y1yLjMrKiZfWRmK0LnGXayMeVkXBnO+YGZLYuVuMzsFDMbaWarQj+3zjESV5DXfWszWxz6/MVm1jLsXEFe9xHjioHrPtOfV9g5T/y6d87F5QMoD9QLPS8BrAKqA08C/UPb+wNPhJ5XB5YCpwJVgLVAQfyaGZuBM0P7PQkMysW4igHNgFuBFzOcayHQGDDgY+DyoOMCigItQs9PAebEQlxh57saeAtYlsvX17H+HQcDj4aeF0i71gL+dwz6uq8LVAg9rwH8HCPXfcS4YuC6z/TnlZ3r/qT/g+S1B/A+fh3slUD5sH+ElaHnA4ABYft/EroICwNbgHNCF+QrQO/ciitsvxsz/AcuD6wIe90dGBF0XBHO8zzQKxbiAooDc0P/0bKVIHI4rv8BxXIynhy4vmLiug9tN2Ab/o+3mLjuM8YV4b1ArvtIcWXnuo/bJqZwZpaAz7ALgHLOuY0Aoa9lQ7tVxP9HTbMeqOicOwjcBnwLbMD/kEflYlyZqRiK8Yh4YyCu8POcDlwJzIiRuIYCTwN7cyKenIgr9DMCGGpmX5nZO2ZWLui4Yuy67wx87ZzbT2xd9+FxhZ/ndIK97jPGddLXfdwnCDMrDrwL3OOc23WsXSNsc2ZWGP8fpS5QAfgGf7eRW3FleooI27I9ZjkH4ko7TyFgPPCCc25d0HGZWR3gT865SdmNJSfjwjflnA3Mc87VA74Ango6rli57s3sQuAJoE/apgi75fp1HyGutO2BXvcZ48rudR/XCSJ0kb8LjHPOpS0+vcnMyofeL49vZwX/l0ilsMPPxv/lVAfAObfW+fu1/wBNcjGuzKwPxZgx3qDjSjMSWO2cey47MeVgXI2B+maWjL/dPs/MZsVAXNvwf9ml/Qd+B6gXA3HVgWCvezM7G/9zud45tza0OfDrPpO40gR23WcSV7au+7hNEKGRDaOA5c65Z8Le+gC4IfT8BnzbXtr2bmZ2qplVAariO8N+BqqbWVq1w9bA8lyMK6LQ7eVuM2sUOuf1xzsmN+IKnetRoBRwz8nGk9NxOeeGO+cqOOcS8J2yq5xzl8ZAXA6YDKTF0gr4Pui4CPi6DzXTfITvF5yXtnPQ131mcYXeC+y6P8bPK3vXfU51osTaI/TDcPhb4yWhR3ugNL5tcHXo6x/CjnkYP3ppJWEjEPAjPJaHzjUZKJ3LcSUD24E9+L+gqoe2JwLLQjG/SGhmfJBx4f+ic6GfV9p5bgk6rgznTCD7o5hy8t/xHGB26FwzgMoxEldg1z3wCPBb2L5LgLJBX/eZxUXA1/2xfl7Zue5VakNERCKK2yYmERHJHiUIERGJSAlCREQiUoIQEZGIlCBERCQiJQgREYlICUIkhphZwaBjEEmjBCFyksxsaFqd/tDrf5jZ3WZ2v5ktMr+OwuCw998L1er/zsx6h23fY2ZDzGwBvjSCSExQghA5eaMIlT0wswJAN2ATvkxLA3w9o/pmdnFo/5udc/XxM4HvNrPSoe3F8DNcGzrn5uZi/CLHVCjoAETyKudcspltM7O6QDnga+AioE3oOfha/FXxpTTuNrNOoe2VQtu3ASn4omwiMUUJQiR7XsMvtnMWMBpfbO9x59yI8J3M7FLgMqCxc25vqKJmkdDb+5xzKbkUr0iWqYlJJHsmAe3wdw6fhB43h+r4Y2YVzawsvsrnr6HkUA1oFFTAIlmlOwiRbHDOHTCzz4AdobuAaWZ2AfCFr9jMHqAnMBW41cy+wVcL/jKomEWyStVcRbIh1Dn9FdDFObc66HhEcpKamEROkplVB9YAM5QcJB7pDkJERCLSHYSIiESkBCEiIhEpQYiISERKECIiEpEShIiIRPT/IpvIeh+KkGgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    " #预测\n",
    "def GM11(x,n):\n",
    "    '''\n",
    "    灰色预测\n",
    "    x：序列，numpy对象\n",
    "    n:需要往后预测的个数\n",
    "    '''\n",
    "    x1 = x.cumsum()#一次累加  \n",
    "    z1 = (x1[:len(x1) - 1] + x1[1:])/2.0#紧邻均值  \n",
    "    z1 = z1.reshape((len(z1),1))  \n",
    "    B = np.append(-z1,np.ones_like(z1),axis=1)  \n",
    "    Y = x[1:].reshape((len(x) - 1,1))\n",
    "    #a为发展系数 b为灰色作用量\n",
    "    [[a],[b]] = np.dot(np.dot(np.linalg.inv(np.dot(B.T, B)), B.T), Y)#计算待估参数  \n",
    "    result = (x[0]-b/a)*np.exp(-a*(n-1))-(x[0]-b/a)*np.exp(-a*(n-2))  #预测方程\n",
    "    S1_2 = x.var()#原序列方差\n",
    "    e = list()#残差序列\n",
    "    for index in range(1,x.shape[0]+1):\n",
    "        predict = (x[0]-b/a)*np.exp(-a*(index-1))-(x[0]-b/a)*np.exp(-a*(index-2))\n",
    "        e.append(x[index-1]-predict)\n",
    "        print(predict)    #预测值\n",
    "    S2_2 = np.array(e).var()#残差方差\n",
    "    C = S2_2/S1_2#后验差比\n",
    "    if C<=0.35:\n",
    "        assess = '后验差比<=0.35，模型精度等级为好'\n",
    "    elif C<=0.5:\n",
    "        assess = '后验差比<=0.5，模型精度等级为合格'\n",
    "    elif C<=0.65:\n",
    "        assess = '后验差比<=0.65，模型精度等级为勉强'\n",
    "    else:\n",
    "        assess = '后验差比>0.65，模型精度等级为不合格'\n",
    "    #预测数据\n",
    "    predict = list()\n",
    "    for index in range(x.shape[0]+1,x.shape[0]+n+1):\n",
    "        predict.append((x[0]-b/a)*np.exp(-a*(index-1))-(x[0]-b/a)*np.exp(-a*(index-2)))\n",
    "        #print((x[0]-b/a)*np.exp(-a*(index-1)))\n",
    "        #print((x[0]-b/a)*np.exp(-a*(index-2)))\n",
    "    predict = np.array(predict)\n",
    "  \n",
    "\n",
    "    return {\n",
    "            'a':{'value':a,'desc':'发展系数'},\n",
    "            'b':{'value':b,'desc':'灰色作用量'},\n",
    "            'predict':{'value':result,'desc':'第%d个预测值'%n},\n",
    "            'C':{'value':C,'desc':assess},\n",
    "            'predict':{'value':predict,'desc':'往后预测%d个的序列'%(n)},\n",
    "            }\n",
    " \n",
    "if __name__ == \"__main__\":\n",
    "    data = np.array([109276.2,120480.4,136576.3,161415.4,185998.9,219028.5,270704,321229.5,347934.9,410354.1,483392.8,537329,588141.2,644380.2,686255.7,743408.3,831381.2,914327.1,984179,1009151])\n",
    "    x = data[0:20]#输入数据\n",
    "    y = data[0:]#需要预测的数据\n",
    "    result = GM11(x,len(y))\n",
    "    predict = result['predict']['value']\n",
    "    predict = np.round(predict,1)\n",
    "    print('真实值:',y)\n",
    "    print('预测值:',predict)\n",
    "    print(result)\n",
    "#作图\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "x1 = np.array([2008,2009,2010,2011,2012,2013,2014,2015,2016,2017,2018,2019])\n",
    "y1 = np.array([7.19,8.26,9.5,9.98,10.48,11,12.5,13.7,14.79,15.66,16.3,16.66])\n",
    "x2 = np.array([2013,2014,2015,2016,2017,2018,2019,2020,2021,2022,2023,2024])\n",
    "y2 = np.array([11.4,12.3,13.3,14.3,15.3,16.5,17.8,19.1,20.6,22.2,23.9,25.7])\n",
    "plt.plot(x1,y1,'r*-',label='true value')   #真实值\n",
    "plt.plot(x2,y2,'b+-',label='predicted value')   #预测值\n",
    "plt.xlabel('year')\n",
    "plt.ylabel('value')\n",
    "plt.legend()\n",
    "plt.plot()\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}